Long range Coulomb interactions in bilayer graphene 
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We report on our studies of interacting electrons in bilayer graphene in a magnetic field. We 
demonstrate that the long range Coulomb interactions between electrons in this material are highly 
important and account for the band asymmetry in recent optical magneto-absorption experiments 
We show that in the unbiased bilayer (where both layers are at the same electrostatic potential), 
the interactions can cause mixing of Landau levels in moderate magnetic fields. For the biased 
bilayer (when the two layers are at different potentials), we demonstrate that the interactions are 
responsible for a change in the total spin of the ground state for half-filled Landau levels in the 
valence band. 



Monolayer graphene is a two-dimensional hexagonal 
crystal of carbon atoms whose gapless, relativistic-like, 
linear low energy dispersion has made it the subject of 
intense study since it was first isolated in 2004 [2]. Bi- 
layer graphene (BLG) the subject of our present 
study consists of a pair of monolayers bound by relatively 
weak dimer bonds formed perpendicular to the mono- 
layer planes. Both the conduction and valence bands 
have low energy structure consisting of two quadratic 
branches separated by the energy associated with the 
dimer bond, 71, and the lower conduction band and up- 
per valence band are degenerate at the K points of the 
Brillouin zone. The existence of chiral charge carriers 
with a Berry's phase of 27r was confirmed in the obser- 
vation of the integer quantum Hall effect [3| where the 
low energy Landau level (LL) spectrum is approximately 
linear in the field with £"„ ~ -iihujc\/n{n — 1) for n > 
where Wc is the cyclotron frequency, and the spectrum 
includes a doubly degenerate LL at zero energy It 
has been predicted theoretically Baa and observed 
experimentally [1, @] , that a gap can be induced in the 
low energy band structure by breaking the symmetry be- 
tween the two layers. Switching of the conduction current 
by sweeping the Fermi energy through the gapped region 
has been observed at low temperatures [lo!| , and this has 
lead to a surge of interest in gapped BLG. 

While the single particle theory of BLG is well known 
0, H, d, 0] , it has been shown that the electron-electron 
interactions are significant in monolayer graphene [ill . 
[l2l |. The Coulomb interaction (CI) has been studied in 
the ungapped bilayer [l3| , while the biasing potential was 
considered in the context of a ferromagnetic transition 
due to short-range interactions in the mean-field approx- 
imation at zero magnetic field [l^, and the absence of 
contribution to the intra-LL cyclotron resonance from the 
electron interactions has been predicted within Hartree- 
Fock theory [l^. However, the effect of the long range 
CI on the ground state of the biased system in a mag- 
netic field has not been investigated, and we address this 
problem in the current Letter. 

We find that the long range nature of the CI makes 
significant changes to the properties of the low energy 
charge carriers for BLG in a magnetic field. The intcrac- 



FIG. 1; The lattice structure of bilayer graphene. The up- 
per (lower) lattice is shown by solid (dashed) lines, (a) The 
top-down view; (b) the side-on view projected along the axis 
between the two arrows in (a). 



tions are significantly stronger for electrons in the lowest 
LL, and this manifests itself in an observable way by lift- 
ing the degeneracy of the cyclotron resonance transitions 
at filling factors ±4 P, [l^ . It also allows the possibility 
of mixing of the LLs which were well-separated in energy 
when the CIs were not considered. By calculating the 
explicit form of the ground state wave function, we show 
how this mixing fundamentally changes the nature of the 
ground state in the biased bilayer, by inducing a finite 
spin polarization for half-filled LLs. 

We model BLG as two sheets, each containing two in- 
equivalent triangular sublattices (labelled A and B) of 
carbon atoms. In the Bernal stacking arrangement, the 
inter-layer bonds consist of dimers formed from atomic 
orbitals associated with the A sublattice in one layer and 
the B sublattice in the other (see Fig. [1]) , and the energy 
associated with this bond is denoted 71 [l3|. We allow 
for a static electric potential U to be applied between 
the layers, so that the upper (lower) layer has potential 
C//2 (— C//2). In a strong magnetic field we can write the 
tight-binding Hamiltonian using a four component single 
valley basis where ^ = ±1 labels the valley as 
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where tt and tt'^ are the operators corresponding to elec- 
tron hops between neighboring atoms (in opposite sublat- 
tices in the same layer). The spectrum is found from 
the quartic polynomial derived from the Schrodinger 
equation associated with the Hamiltonian in Eq. ([1]), 
where n E {0, 1,2,.. .}. Additionally, we denote the band 
of a particular LL by placing a '+' (' — ') after the level 
index for the conduction (valence) band. The wave func- 
tions associated with the LLs from the two low energy 
bands are given by 
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where the functions are the magnetic oscillator func- 
tions in the Landau gauge, and the coefficients a, 6, c, 
and d are defined so that the overall wave function is 
normalized to unity. There are also levels with n = 0±, 
which have wave functions V'o+ = £''^^('(^'0,0,0,0), and 
^ifj ^ f^ky ((^^^0,0,0) with £ = ±5] and and 
are defined as higher LLs above with the appropriate sub- 
stitutions for n and ^. When ?7 = 0, these four states are 
degenerate yielding the eight-fold degeneracy (including 
the factor of 2 for spin) seen in the integer quantum Hall 
effect in BLG 4]. We include the fermionic properties of 
the electrons by constructing Slater determinants for the 
non-interacting many body basis wave functions. 

To include the effects of the long-range CI we consider 
the Hamiltonian 



2 ^ e|r,; 
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where the vectors rj.j label the positions of the electrons, 
and e = Aireon is the dielectric constant of graphene. For 
graphene mounted on an Si02 substrate, k = 2.5 [l^. 
Our analysis is conducted by employing the exact diag- 
onalization scheme [l9l | in which we calculate the linear 
combination of non-interacting basis states which gives 
the ground state of the Hamiltonian H = Ho + Hcoui- 
This method entails dividing the infinite sheet into rect- 
angular cells of dimension x Ly [20i\ and applying pe- 
riodic boundary conditions to the wavefunctions at the 
edges of each cell. The matrix elements of the interaction 
over the single particle states are evaluated exactly, and 
the interaction between the cells is taken into account by 
adding the Madelung energy of a charged lattice [2l[ . 

The single particle states included in the Hilbert space 
from which the non-interacting many body basis is con- 
structed are as follows. There are four relevant quantum 
numbers: The LL index n, the valley ^, the spin and 
the momentum k = Trm/Lx- The values of the momen- 
tum are fixed when the boundary conditions are applied 
to the cell, and are labelled byO<TO<M— 1 with 
M — LxLy/{2TTX'^). The LLs selected are governed by 
the details of the system we wish to model, and M is 
set by computational restraints. Our model includes all 
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TABLE I: Energy shift per electron due to the Coulomb in- 
teraction for integer filling factors in the (a) 0± LL and (b) 
\n\ > 1 LLs, for U = 0. Energy units are e^/(eAs), the 
number of momentum states M = 3, and the magnetic field 
B = 3T. 



inter-electron screening effects since we calculate the ex- 
act matrix elements of the full Coulomb interaction, and 
it is well known that filled Landau levels with energy sig- 
nificantly below the Fermi level do not make additional 
contributions. 

In order to reduce the size of the many body sys- 
tem (and so improve the calculation speed), we see that 
the Hamiltonian conserves the total momentum fj, = 
'l2iLi "^i mod M . Therefore, we can perform a seper- 
ate diagonalization for each value of /i, and reduce the 
basis size to approximately the 1/Afth part. We define 
S = Si to be the total spin of the many electron sys- 
tem, and since there is no spin-dependent term in the 
Hamiltonian, Sz (the projection of S on the z axis) is a 
good quantum number. Therefore we fix Sz to its mini- 
mum value whilst still being able to recover all eigenstates 
of S'^ 12 1, further reducing the many body basis size. 

We consider two cases: Firstly we demonstrate the 



strength of the interaction by calculating the shift in the 
energy of each LL due to interactions for [/ = (an unbi- 
ased bilayer) and apply the results to recent experimental 
data. Then we examine the system where the filling fac- 
tor is negative, the inter-layer potential sizeable, and the 
magnetic field strong. In this case, we observe changes 
in the total spin of the ground state as a function of U 
and the magnetic field strength B. 

We model the unbiased bilayer near half-filling by tak- 
ing a single particle Hilbert space consisting of electrons 
in the 0-1- and 0— LLs with all possible spin and val- 
ley states at [/ = 0. Each integer value of the filling 
factor V is simulated by taking the number of electrons 
N = {v + A)M, and we have M = 3. Tablelja) shows the 
results of diagonalizing the resulting many body Hamilto- 
nian and evaluating the change in energy from the non- 
interacting ground state for integer filling factors. We 
notice that the energy shift per electron reduces slightly 
as the LL fills. 

In Table mb), we show the energy shift due to the CI 
for electrons in higher LLs (i.e. for levels with \n\ > 1). 



3 



100 

I 60 

S 20 

c 80 
o 

1 ^° 

2 40 

20- 
0, 



(a) 71 


0.40eV 




. 




- 


- (b) 71 = 


0.35eY. 




1/ = -4 














— s = 


1.00- 






0.95 




€> — 


0.90 



100 



10 15 
Magnetic field [T] 



20 5 




Electrons in filled Landau levels 



10 15 
Magnetic field [T] 



FIG. 2: Electron-hole asymmetry in the inter-LL optical 
transition energy. The experimental data (represented as 
points) are taken from Ref. Q Fig. 2. In (a) and (b), 
V — «F/(10^ms~'^); in (c) we take 71 = 0.4eV, vf — 
0.95 X lO^'ms-^ M = 6 throughout. 



We have taken a single particle Hilbert space consisting of 
all spin and valley states within one LL. The filling factor 
Vn within LL n can range between (corresponding to 
an empty level) and 4 (a filled level), so that Vn — 4 
and = describe the same overall filling factor. 

The number of electrons is set by A'^ = VnM, and in 
order to allow direct comparison with the lowest LL we 
restrict ourselves to M = 3. The energy associated with 
the interaction of electrons is very similar in each of the 
higher LLs, and that the interaction energy per particle 
is slightly reduced as the LL is filled. We have verified 
that the results are identical in the valence band. 

Together, Tables|Tl[a,b) show that the effect of the long 
range CIs is considerable, and that for this value of the 
magnetic field {B = 3T) the shift in the higher LLs is only 
about two-thirds that of the lowest LL. This difference in 
the shift will reveal itself in the infra-red absorption spec- 
trum of bilayer graphene, since the energy of the optical 
transitions depends entirely on the direct energy spacing 
between levels. At U ~ and = —4, the lowest energy 
transition is from the 1— level to the 0± level, while at 
1/ — +4 the lowest energy transition is from 0± to H-. 
Therefore, if the 0± is shifted with respect to the two 
\n\ = 1 levels, the degeneracy of these two transitions 
predicted in the single particle theory [l^ will be lifted. 
In Figure [5] we show the predictions of our theory in com- 
parison to recent experimental data [if. Panes (a) and 

(b) show comparison of the valence band transitions {i.e. 
for 1/ = —4) for two values of 71 and three values of vp- 
The energy of the excited state is calculated as the sum 
of the full interacting energy of — 1 electrons in the 
1— level and 1 electron in the 0± level (where we assume 
no LL mixing because the magnetic field is strong). The 
transition energy is the difference between this and the 
energy of N interacting electrons in the 1— level. Pane 

(c) shows the comparison of theoretical and experimen- 
tal data in both bands for the best values of parameters. 
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FIG. 3: (a) The energy shift per electron of filled LLs. (b) 
The absolute energy per electron of filled LLs showing the 
crossing between the n = 0± degenerate level and the higher 
LLs in the valence band. In both plots U = 0, and M = 5. 



The correspondence to our theory is clear. 

In Fig. [3] we show the energy shift and absolute energy 
of filled LLs as a function of the magnetic field. The 



strength of the interaction scales with e^/(eAB) cx vi? 
with a roughly constant coefficient, while the LL spacing 
goes as huc oc B, so for lower values of the field, the 
n = 0± level crosses the 2— and 1— levels as shown in 
Fig. [3jb). The data shown here were calculated with 
K = 2.5, modelHng graphene deposited on an Si02 
substrate. For suspended graphene (where k « 1), it is 
conceivable that the effect of the interaction would be 
even stronger. Additionally, the effect of the inter-layer 
potential is to bring together the valence band LLs with 
low index Q , so it is plausible that the interactions will 
cause significant mixing between these levels. 

Now we turn our attention to the system with an inter- 
layer potential applied, so that U ^ 0. With a finite gap 
between the 0-1- and 0— levels and non-zero filling fac- 
tor, it is possible to consider the negatively-doped system 
by taking only those single particle states which are in 
the valence band. Therefore we select the single particle 
states which form the Slater determinants describing the 
non-interacting basis states by taking all spin and val- 
ley states of the 0— and 1— LLs. We have M — 2 and 
the number of electrons is related to the filling factor by 
N — {ly + 8)M. Diagonalizing this system for half-filled 
LLs (so for v = —6 and ^ = —2), and calculating the 
expectation value of the total spin operator S'^ over the 
resulting ground state as a function of the magnetic field 
and the inter-layer bias gives the data shown in the plots 
in Fig. |4l^a,b). We have superimposed lines which repre- 
sent the energy at which the single particle energy levels 
cross, as labelled in the legend. We have also calculated 
the expectation of the total valley quantum number for 
each of these systems and find that it is constant with 
a value equal to 2. This is expected for half- filled LLs 
because of the lifting of the valley degeneracy in the sin- 
gle particle theory by the inter-layer potential @, 0] ■ The 
plots show that there is an abrupt change in the total 
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dated with the 0± LL, and LL mixing between the 0± 
and 1— levels is induced for moderate magnetic fields. If 
an inter-layer bias is applied to split the valence and con- 
duction bands, the electron-electron interactions precip- 
itate a transition in the total spin of the ground state of 
half-filled LLs for certain ranges of parameters. Various 
experimental techniques to measure ground state spin 
polarizations in the quantum Hall effect regime were re- 
ported earlier 22| and might prove to be useful here as 
well. This effect will have fundamental implications for 
the design of devices made from this material. 

We acknowledge support from the Canada Research 
Chairs Program, and the NSERC Discovery Grant. 



FIG. 4: The total spin of the ground state of the (a) v — —2 
and (b) ty ^ -6 systems. M = 2 and N = + 8)M. The 
lines show the crossing points of the single particle states. The 
graining is due to the finite interval between data points, (c) 
The occupancy of the single electron states in the interacting 
many body ground state for each region of the plot in (a). 
The z-projection of the total spin is fixed at zero as described 
in the text. 



spin of the ground state, and a range of parameters where 
there is a non-zero polarization of the spin. This transi- 
tion is not directly related to the crossings of the single 
particle states, since the position and slope of the transi- 
tions do not match the corresponding lines superimposed 
on the plots. This effect is therefore purely due to the 
CI, and in particular to the exchange contribution which 
acts to minimize the energy of spin-polarized many body 
states. 

Figure IH^c) shows the occupation of the single particle 
levels in the interacting many body ground state of the 
V = —2 system. For simplicity, we display only the LL 
index of the states. The actual ground state is a coher- 
ent combination of several non-interacting basis states, 
where the combination of LL indices is consistent but 
different arrangements of momentum and valley states 
each come with identical prefactors in the linear combi- 
nation. In the lower-right region of the parameter space, 
the electrons occupy as many of the 1— states as possi- 
ble. Moving toward the upper-left region, the 0— levels 
become successively more populated. The absence of spin 
polarization in regions 1 and 3 is caused by the pairing of 
electrons in the same valley. In region two, where there 
are six electrons per LL, this pairing is incomplete and 
the spin polarization finite. The pattern of filling in the 
two regions of the — —6 system is identical. 

In conclusion, we have shown that the long range CI 
between electrons plays an important role near the Dirac 
point in BLG. In the unbiased case, the interactions will 
cause a change in the cyclotron resonance energies asso- 
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